{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# C/C++ Library with Python \n",
    "\n",
    "* Change the current working directory into `./demo`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "J:\\SEU\\SEE\\PySEE\\home\\notebook\\demo\n"
     ]
    }
   ],
   "source": [
    "%cd demo"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "'J:\\\\SEU\\\\SEE\\\\PySEE\\\\home\\\\notebook'"
      ]
     },
     "execution_count": 1,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "%pwd"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "\n",
    "The advantage of Python is that it is **flexible and easy** to program. The time it takes to setup a new calulation is therefore short. \n",
    "\n",
    "But for certain types of calculations Python (and any other interpreted language) can be **very slow**.\n",
    "\n",
    "Such calculations may be implemented in a compiled language such as C or Fortran."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 1 ctypes - access the C library\n",
    "\n",
    "**ctypes** is a foreign function library for Python. It provides C compatible data types, and allows calling functions in DLLs or shared libraries. It can be used to wrap these libraries in pure Python.\n",
    "\n",
    "http://docs.python.org/3/library/ctypes.html\n",
    "\n",
    "We need to load the library and set properties such as the functions return and argument types using the **ctypes** package :\n",
    "\n",
    "1. **Loads the library** \n",
    "\n",
    "  **ctypes** exports the **cdll**, and on Windows **windll** objects, for loading dynamic link libraries.\n",
    "\n",
    "   * **cdll.LoadLibrary(name)** : loads the library  which export functions using standard `__cdecl` calling convention\n",
    "   * **windll.LoadLibrary(name)** : loads the library with `__stdcall` calling convention for the function \n",
    "\n",
    "\n",
    "2. **Specifying the required `argument` types:`argtypes`**  \n",
    "  \n",
    "  * It specify the required argument types of functions exported from DLLs by setting the argtypes attribute\n",
    "\n",
    "3. **Return types: `restype`**\n",
    "\n",
    "  * Return typescan be specified by setting the restype attribute of the function object.\n",
    "\n",
    "4. **Function prototypes**\n",
    " \n",
    " * **CFUNCTYPE** \n",
    " \n",
    "   The CFUNCTYPE factory function creates types for callback functions using the normal cdecl calling convention\n",
    "   \n",
    "```python\n",
    "ctypes.CFUNCTYPE(restype, *argtypes, use_errno=False, use_last_error=False)\n",
    "``` \n",
    " * **WINFUNCTYPE**\n",
    "  \n",
    "  Windows only: The returned function prototype creates functions that use the `__stdcall` calling convention\n",
    "```python\n",
    "ctypes.WINFUNCTYPE(restype, *argtypes, use_errno=False, use_last_error=False)\n",
    "```"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 2 Call the Shared Library from Python\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 2.1  C compatible data types\n",
    "\n",
    "#### 2.1.1 Fundamental data types\n",
    "\n",
    "ctypes defines a number of primitive C compatible data types:\n",
    "\n",
    "https://docs.python.org/3/library/ctypes.html#fundamental-data-types\n",
    "\n",
    "\n",
    "\n",
    "<table class=\"docutils\" border=\"1\">\n",
    "<colgroup>\n",
    "<col width=\"24%\">\n",
    "<col width=\"46%\">\n",
    "<col width=\"30%\">\n",
    "</colgroup>\n",
    "<thead valign=\"bottom\">\n",
    "<tr class=\"row-odd\"><th class=\"head\">ctypes type</th>\n",
    "<th class=\"head\">C type</th>\n",
    "<th class=\"head\">Python type</th>\n",
    "</tr>\n",
    "</thead>\n",
    "<tbody valign=\"top\">\n",
    "<tr class=\"row-even\"><td><a title=\"ctypes.c_bool\" class=\"reference internal\" href=\"#ctypes.c_bool\"><code class=\"xref py py-class docutils literal notranslate\"><span class=\"pre\">c_bool</span></code></a></td>\n",
    "<td><code class=\"xref c c-type docutils literal notranslate\"><span class=\"pre\">_Bool</span></code></td>\n",
    "<td>bool (1)</td>\n",
    "</tr>\n",
    "<tr class=\"row-odd\"><td><a title=\"ctypes.c_char\" class=\"reference internal\" href=\"#ctypes.c_char\"><code class=\"xref py py-class docutils literal notranslate\"><span class=\"pre\">c_char</span></code></a></td>\n",
    "<td><code class=\"xref c c-type docutils literal notranslate\"><span class=\"pre\">char</span></code></td>\n",
    "<td>1-character bytes object</td>\n",
    "</tr>\n",
    "<tr class=\"row-even\"><td><a title=\"ctypes.c_wchar\" class=\"reference internal\" href=\"#ctypes.c_wchar\"><code class=\"xref py py-class docutils literal notranslate\"><span class=\"pre\">c_wchar</span></code></a></td>\n",
    "<td><code class=\"xref c c-type docutils literal notranslate\"><span class=\"pre\">wchar_t</span></code></td>\n",
    "<td>1-character string</td>\n",
    "</tr>\n",
    "<tr class=\"row-odd\"><td><a title=\"ctypes.c_byte\" class=\"reference internal\" href=\"#ctypes.c_byte\"><code class=\"xref py py-class docutils literal notranslate\"><span class=\"pre\">c_byte</span></code></a></td>\n",
    "<td><code class=\"xref c c-type docutils literal notranslate\"><span class=\"pre\">char</span></code></td>\n",
    "<td>int</td>\n",
    "</tr>\n",
    "<tr class=\"row-even\"><td><a title=\"ctypes.c_ubyte\" class=\"reference internal\" href=\"#ctypes.c_ubyte\"><code class=\"xref py py-class docutils literal notranslate\"><span class=\"pre\">c_ubyte</span></code></a></td>\n",
    "<td><code class=\"xref c c-type docutils literal notranslate\"><span class=\"pre\">unsigned</span> <span class=\"pre\">char</span></code></td>\n",
    "<td>int</td>\n",
    "</tr>\n",
    "<tr class=\"row-odd\"><td><a title=\"ctypes.c_short\" class=\"reference internal\" href=\"#ctypes.c_short\"><code class=\"xref py py-class docutils literal notranslate\"><span class=\"pre\">c_short</span></code></a></td>\n",
    "<td><code class=\"xref c c-type docutils literal notranslate\"><span class=\"pre\">short</span></code></td>\n",
    "<td>int</td>\n",
    "</tr>\n",
    "<tr class=\"row-even\"><td><a title=\"ctypes.c_ushort\" class=\"reference internal\" href=\"#ctypes.c_ushort\"><code class=\"xref py py-class docutils literal notranslate\"><span class=\"pre\">c_ushort</span></code></a></td>\n",
    "<td><code class=\"xref c c-type docutils literal notranslate\"><span class=\"pre\">unsigned</span> <span class=\"pre\">short</span></code></td>\n",
    "<td>int</td>\n",
    "</tr>\n",
    "<tr class=\"row-odd\"><td><a title=\"ctypes.c_int\" class=\"reference internal\" href=\"#ctypes.c_int\"><code class=\"xref py py-class docutils literal notranslate\"><span class=\"pre\">c_int</span></code></a></td>\n",
    "<td><code class=\"xref c c-type docutils literal notranslate\"><span class=\"pre\">int</span></code></td>\n",
    "<td>int</td>\n",
    "</tr>\n",
    "<tr class=\"row-even\"><td><a title=\"ctypes.c_uint\" class=\"reference internal\" href=\"#ctypes.c_uint\"><code class=\"xref py py-class docutils literal notranslate\"><span class=\"pre\">c_uint</span></code></a></td>\n",
    "<td><code class=\"xref c c-type docutils literal notranslate\"><span class=\"pre\">unsigned</span> <span class=\"pre\">int</span></code></td>\n",
    "<td>int</td>\n",
    "</tr>\n",
    "<tr class=\"row-odd\"><td><a title=\"ctypes.c_long\" class=\"reference internal\" href=\"#ctypes.c_long\"><code class=\"xref py py-class docutils literal notranslate\"><span class=\"pre\">c_long</span></code></a></td>\n",
    "<td><code class=\"xref c c-type docutils literal notranslate\"><span class=\"pre\">long</span></code></td>\n",
    "<td>int</td>\n",
    "</tr>\n",
    "<tr class=\"row-even\"><td><a title=\"ctypes.c_ulong\" class=\"reference internal\" href=\"#ctypes.c_ulong\"><code class=\"xref py py-class docutils literal notranslate\"><span class=\"pre\">c_ulong</span></code></a></td>\n",
    "<td><code class=\"xref c c-type docutils literal notranslate\"><span class=\"pre\">unsigned</span> <span class=\"pre\">long</span></code></td>\n",
    "<td>int</td>\n",
    "</tr>\n",
    "<tr class=\"row-odd\"><td><a title=\"ctypes.c_longlong\" class=\"reference internal\" href=\"#ctypes.c_longlong\"><code class=\"xref py py-class docutils literal notranslate\"><span class=\"pre\">c_longlong</span></code></a></td>\n",
    "<td><code class=\"xref c c-type docutils literal notranslate\"><span class=\"pre\">__int64</span></code> or <code class=\"xref c c-type docutils literal notranslate\"><span class=\"pre\">long</span> <span class=\"pre\">long</span></code></td>\n",
    "<td>int</td>\n",
    "</tr>\n",
    "<tr class=\"row-even\"><td><a title=\"ctypes.c_ulonglong\" class=\"reference internal\" href=\"#ctypes.c_ulonglong\"><code class=\"xref py py-class docutils literal notranslate\"><span class=\"pre\">c_ulonglong</span></code></a></td>\n",
    "<td><code class=\"xref c c-type docutils literal notranslate\"><span class=\"pre\">unsigned</span> <span class=\"pre\">__int64</span></code> or\n",
    "<code class=\"xref c c-type docutils literal notranslate\"><span class=\"pre\">unsigned</span> <span class=\"pre\">long</span> <span class=\"pre\">long</span></code></td>\n",
    "<td>int</td>\n",
    "</tr>\n",
    "<tr class=\"row-odd\"><td><a title=\"ctypes.c_size_t\" class=\"reference internal\" href=\"#ctypes.c_size_t\"><code class=\"xref py py-class docutils literal notranslate\"><span class=\"pre\">c_size_t</span></code></a></td>\n",
    "<td><code class=\"xref c c-type docutils literal notranslate\"><span class=\"pre\">size_t</span></code></td>\n",
    "<td>int</td>\n",
    "</tr>\n",
    "<tr class=\"row-even\"><td><a title=\"ctypes.c_ssize_t\" class=\"reference internal\" href=\"#ctypes.c_ssize_t\"><code class=\"xref py py-class docutils literal notranslate\"><span class=\"pre\">c_ssize_t</span></code></a></td>\n",
    "<td><code class=\"xref c c-type docutils literal notranslate\"><span class=\"pre\">ssize_t</span></code> or\n",
    "<code class=\"xref c c-type docutils literal notranslate\"><span class=\"pre\">Py_ssize_t</span></code></td>\n",
    "<td>int</td>\n",
    "</tr>\n",
    "<tr class=\"row-odd\"><td><a title=\"ctypes.c_float\" class=\"reference internal\" href=\"#ctypes.c_float\"><code class=\"xref py py-class docutils literal notranslate\"><span class=\"pre\">c_float</span></code></a></td>\n",
    "<td><code class=\"xref c c-type docutils literal notranslate\"><span class=\"pre\">float</span></code></td>\n",
    "<td>float</td>\n",
    "</tr>\n",
    "<tr class=\"row-even\"><td><a title=\"ctypes.c_double\" class=\"reference internal\" href=\"#ctypes.c_double\"><code class=\"xref py py-class docutils literal notranslate\"><span class=\"pre\">c_double</span></code></a></td>\n",
    "<td><code class=\"xref c c-type docutils literal notranslate\"><span class=\"pre\">double</span></code></td>\n",
    "<td>float</td>\n",
    "</tr>\n",
    "<tr class=\"row-odd\"><td><a title=\"ctypes.c_longdouble\" class=\"reference internal\" href=\"#ctypes.c_longdouble\"><code class=\"xref py py-class docutils literal notranslate\"><span class=\"pre\">c_longdouble</span></code></a></td>\n",
    "<td><code class=\"xref c c-type docutils literal notranslate\"><span class=\"pre\">long</span> <span class=\"pre\">double</span></code></td>\n",
    "<td>float</td>\n",
    "</tr>\n",
    "<tr class=\"row-even\"><td><a title=\"ctypes.c_char_p\" class=\"reference internal\" href=\"#ctypes.c_char_p\"><code class=\"xref py py-class docutils literal notranslate\"><span class=\"pre\">c_char_p</span></code></a></td>\n",
    "<td><code class=\"xref c c-type docutils literal notranslate\"><span class=\"pre\">char</span> <span class=\"pre\">*</span></code> (NUL terminated)</td>\n",
    "<td>bytes object or <code class=\"docutils literal notranslate\"><span class=\"pre\">None</span></code></td>\n",
    "</tr>\n",
    "<tr class=\"row-odd\"><td><a title=\"ctypes.c_wchar_p\" class=\"reference internal\" href=\"#ctypes.c_wchar_p\"><code class=\"xref py py-class docutils literal notranslate\"><span class=\"pre\">c_wchar_p</span></code></a></td>\n",
    "<td><code class=\"xref c c-type docutils literal notranslate\"><span class=\"pre\">wchar_t</span> <span class=\"pre\">*</span></code> (NUL terminated)</td>\n",
    "<td>string or <code class=\"docutils literal notranslate\"><span class=\"pre\">None</span></code></td>\n",
    "</tr>\n",
    "<tr class=\"row-even\"><td><a title=\"ctypes.c_void_p\" class=\"reference internal\" href=\"#ctypes.c_void_p\"><code class=\"xref py py-class docutils literal notranslate\"><span class=\"pre\">c_void_p</span></code></a></td>\n",
    "<td><code class=\"xref c c-type docutils literal notranslate\"><span class=\"pre\">void</span> <span class=\"pre\">*</span></code></td>\n",
    "<td>int or <code class=\"docutils literal notranslate\"><span class=\"pre\">None</span></code></td>\n",
    "</tr>\n",
    "</tbody>\n",
    "</table>\n",
    " \n",
    "All these types can be created by calling them with an optional initializer of the correct type and value:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from ctypes import *\n",
    "v=c_int(10)\n",
    "# the instance of c_int\n",
    "print(v)\n",
    "# get the value of the instance\n",
    "print(v.value)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "####  2.1.2  Using ctypes to call shared lib\n",
    "\n",
    "We  use `ctypes` call the shared `libfuns.dll`\n",
    "\n",
    "* `cdll` loads libraries which export functions using the `__cdecl` calling convention \n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "120\n"
     ]
    }
   ],
   "source": [
    "from ctypes import *\n",
    "flib = cdll.LoadLibrary(\"./bin/libmaths.dll\")\n",
    "## 1 argtypes,restype\n",
    "flib.factorial.restype=c_int\n",
    "flib.factorial.argtypes=[c_int]\n",
    "n=5\n",
    "print(flib.factorial(n))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "120\n"
     ]
    }
   ],
   "source": [
    "from ctypes import *\n",
    "\n",
    "flib = cdll.LoadLibrary(\"./bin/libmaths.dll\")\n",
    "\n",
    "## 2 Function prototypes\n",
    "prototype = CFUNCTYPE(c_int,c_int)\n",
    "f = prototype((\"factorial\", flib),)\n",
    "print(f(n))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Wrapper the Shared Library to the Python API** \n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [],
   "source": [
    "from ctypes import *\n",
    "\n",
    "flib = cdll.LoadLibrary(\"./bin/libmaths.dll\")\n",
    "prototype = CFUNCTYPE(c_int,c_int)\n",
    "\n",
    "# ---(p,t) ----------------\n",
    "def factorial(n):\n",
    "    f = prototype((\"factorial\", flib),)\n",
    "    result = f(n)\n",
    "    return result\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "120"
      ]
     },
     "execution_count": 10,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "f= factorial(5)\n",
    "f"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 2.2 One-dimensional array \n",
    "\n",
    "[Arrays](https://docs.python.org/3/library/ctypes.html#arrays) are sequences, containing a fixed number of instances of the same type.\n",
    "\n",
    "**create array types**\n",
    "\n",
    "The recommended way to `create array types` is by\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "<class '__main__.c_double_Array_10'>\n",
      "0.0\n",
      "[0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0]\n"
     ]
    }
   ],
   "source": [
    "# c_double*10\n",
    "narray=(c_double*10)()  \n",
    "print(type(narray))\n",
    "\n",
    "\n",
    "print(narray[2])\n",
    "print(list(narray))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 27,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "<__main__.c_double_Array_3 object at 0x00000252A4892440>\n",
      "2.0\n",
      "[1.0, 2.0, 3.0]\n"
     ]
    }
   ],
   "source": [
    "# c_double*len(list1)\n",
    "list1=[1,2,3]\n",
    "narray=(c_double*len(list1))(*list1)\n",
    "print(narray)\n",
    "\n",
    "\n",
    "print(narray[1])\n",
    "print(list(narray))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "```\n",
    "double  mean(double data[], int size);\n",
    "prototype = CFUNCTYPE(c_double,POINTER(c_double),c_int)\n",
    "```\n",
    "* `double data[] > POINTER(c_double)`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 29,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "3.75\n"
     ]
    }
   ],
   "source": [
    "#  one-dimensional array\n",
    "from ctypes import *\n",
    "\n",
    "flib = cdll.LoadLibrary(\"./bin/libmaths.dll\")\n",
    "prototype = CFUNCTYPE(c_double,POINTER(c_double),c_int)\n",
    "\n",
    "# mean\n",
    "def fmean(values):\n",
    "    size=len(values)\n",
    "    v=(c_double*size)(*values)\n",
    "    f = prototype((\"mean\", flib),)\n",
    "    result = f(v,size)\n",
    "    return result\n",
    "\n",
    "list1=[1,2,3,9]\n",
    "\n",
    "print(fmean(list1))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 2.3 Function "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### 2.3.1 Bisection Methon Shared Library in C\n",
    "\n",
    "* using `__cdecl`calling convention\n",
    "\n",
    "$y=f(x)$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 32,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Overwriting ./include/roots.h\n"
     ]
    }
   ],
   "source": [
    "%%file ./include/roots.h\n",
    "#ifndef ROOTS_H\n",
    "#define ROOTS_H\n",
    "\n",
    "#include <math.h>\n",
    "// function definition\n",
    "typedef double (*fun)(double);\n",
    "\n",
    "int rtbis(fun func,double y,double x1, double x2, double xacc,double *rtb);\n",
    "/*\n",
    "\tThe program uses the bisection method to solve the equation\n",
    "\t\tf(x)-y = 0.\n",
    "\tThe solution is to be in [x1,x2] and it is assumed that\n",
    "\t\t(f(x1)-y)*(f(x2)-y) <= 0.\n",
    "\tThe solution is returned in rtb, and it is to be in error by at most xacc.\n",
    "\t\n",
    "\treturn value is an error indicator.\n",
    "\t  If =0, the solution has been computed satisfactorily.\n",
    "\t  If =1, (f(x1)-y)*(f(x2)-y) was greater than 0, contrary to assumption \n",
    "      If =2, exceeded the maximum number of iteration \n",
    "*/\n",
    "#endif\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### 2.3.2The ways to pass arguments to the function\n",
    "\n",
    "\n",
    "Following are the two ways to **pass arguments** to the function:\n",
    "\n",
    "* Pass by value\n",
    "* Pass by reference\n",
    "\n",
    "##### Pass-by-Value\n",
    "\n",
    "In pass-by-value, a **copy** of argument is created and passed into the function. The invoked function works on the \"clone\", and **cannot modify the original copy**.\n",
    "\n",
    "In C, fundamental types (such as int and double) are passed by value. That is, you cannot modify caller's value inside the function\n",
    "\n",
    "#####  Pass-by-Reference\n",
    "\n",
    "In pass-by-reference, a **reference** of the caller's variable is passed into the function. In other words, the invoked function works on the **same** data. \n",
    "\n",
    "* If the invoked function modifies the parameter, the same caller's copy will be modified as well.\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 33,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Overwriting ./src/roots.c\n"
     ]
    }
   ],
   "source": [
    "%%file ./src/roots.c\n",
    "/*\n",
    "  Numerical Recipes http://numerical.recipes\n",
    "*/ \n",
    "#include <math.h>\n",
    "#include \"roots.h\"\n",
    "\n",
    "int rtbis(fun func,double y,double x1, double x2, double xacc,double *rtb)\n",
    "{\n",
    " /*\n",
    "\tThe program uses the bisection method to solve the equation\n",
    "\t\tf(x)-y = 0.\n",
    "\tThe solution is to be in [x1,x2] and it is assumed that\n",
    "\t\t(f(x1)-y)*(f(x2)-y) <= 0.\n",
    "\tThe solution is returned in rtb, and it is to be in error by at most xacc.\n",
    "\t\n",
    "\treturn value is an error indicator.\n",
    "\t  If =0, the solution has been computed satisfactorily.\n",
    "\t  If =1, (f(x1)-y)*(f(x2)-y) was greater than 0, contrary to assumption \n",
    "      If =2, exceeded the maximum number of iteration \n",
    "*/\n",
    "\tconst int IMAX=100; // the maximum number of iteration\n",
    "    float dx,f,fmid,xmid;\n",
    "\n",
    "\tf=(*func)(x1)-y;\n",
    "\tfmid=(*func)(x2)-y;\n",
    "\tif (f*fmid >= 0.0) // endpoints do not straddle y=0\n",
    "       return 1; \n",
    "    // init the root value: rtb\n",
    "\t*rtb = f < 0.0 ? (dx=x2-x1,x1) : (dx=x1-x2,x2);\n",
    "\tfor (int i=1;i<=IMAX;i++) {\n",
    "\t\tfmid=(*func)(xmid=(*rtb)+(dx *= 0.5))-y;\n",
    "\t\tif (fmid <= 0.0) *rtb=xmid;\n",
    "\t\tif (fabs(dx) < xacc || fmid == 0.0) \n",
    "           return 0;\n",
    "  \t}\n",
    "    // Exceeded the maximum number of iteration\n",
    "    return 2;\n",
    "}"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### 2.3.3 Building the Library"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 34,
   "metadata": {},
   "outputs": [],
   "source": [
    "!gcc -c -O3 -Wall -fPIC -o ./obj/roots.o  ./src/roots.c -I./include\n",
    "!gcc -shared -o ./bin/libroots.dll  ./obj/roots.o"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### 2.3.4 Caller in C\n",
    "\n",
    "$f(x)=x^2$\n",
    "\n",
    "$4=f(x)$\n",
    "\n",
    "```\n",
    "int rtbis(fun func,double y,double x1, double x2, double xacc,double *rtb)\n",
    "\n",
    "double xl, xu, epsilon, root;\n",
    "ier=rtbis(fcn,y,xl, xu, epsilon, &root);\n",
    "```\n",
    "*  &root > double *rtb"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 35,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Overwriting ./src/rtbisApp.c\n"
     ]
    }
   ],
   "source": [
    "%%file ./src/rtbisApp.c\n",
    "#include <stdio.h>\n",
    "#include <string.h>\n",
    "#include <stdlib.h>\n",
    "#include \"roots.h\"\n",
    "\n",
    "double fcn(double x)\n",
    "{\n",
    "\tdouble result;\n",
    "    result = x * x;\n",
    "\treturn result;\n",
    "}\n",
    "\n",
    "int main()\n",
    "{\n",
    "\tdouble xl, xu, epsilon, root;\n",
    "\tint ier;\n",
    "\txl = 0.1;\n",
    "\txu = 3.2;\n",
    "\tepsilon = 0.001;\n",
    "    double y=4;\n",
    "\t// Calculate root\n",
    "\tier=rtbis(fcn,y,xl, xu, epsilon, &root);\n",
    "\t// Print answers\n",
    "\tprintf(\"root = %14.7e  ier = %1d\", root, ier);\n",
    "\treturn 0;\n",
    "}\t\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 36,
   "metadata": {},
   "outputs": [],
   "source": [
    "!gcc -o  ./bin/rtbisApp ./src/rtbisApp.c -L./bin/ -lroots -I./include"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 37,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "root = 1.9996581e+000  ier = 0\n"
     ]
    }
   ],
   "source": [
    "!.\\bin\\rtbisApp "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### 2.3.5  Python API \n",
    "\n",
    "\n",
    "**CDLL**\n",
    "\n",
    "* Instances of this class represent loaded shared libraries. Functions in these libraries use the standard `C calling convention`, and are assumed to return int.\n",
    "\n",
    "```\n",
    "int rtbis(fun func,double y,double x1, double x2, double xacc,double *rtb)\n",
    "\n",
    "rtbis.argtypes = [fun,c_double, c_double,c_double,c_double,POINTER(c_double)]\n",
    "rtbis.restype = c_int\n",
    "```\n",
    "\n",
    "* POINTER(c_double) > double *rtb"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 38,
   "metadata": {},
   "outputs": [],
   "source": [
    "from ctypes import * \n",
    "\n",
    "rtapi = CDLL('.\\\\bin\\\\libroots.dll')\n",
    "fun = CFUNCTYPE(c_double, c_double)\n",
    "\n",
    "def biroot(fn,y,xl,xu,eps): \n",
    "    rtbis = rtapi.rtbis\n",
    "    rtbis.argtypes = [fun,c_double, c_double,c_double,c_double,POINTER(c_double)]\n",
    "    rtbis.restype = c_int\n",
    "\n",
    "    rtb=c_double()\n",
    "    ier=rtbis(fn,y,xl,xu,eps,byref(rtb))\n",
    "    return rtb.value,ier"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 39,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "1.999902367591858 0\n"
     ]
    }
   ],
   "source": [
    "def fn(x):        \n",
    "    return x**2\n",
    "\n",
    "_fn = fun(fn)\n",
    "xl = c_double(0.2) \n",
    "xu = c_double(3)\n",
    "eps=0.001\n",
    "y=4\n",
    "\n",
    "rtb,ier=biroot(_fn,y,xl,xu,eps)\n",
    "print(rtb,ier)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "###  2.4 Call the CoolProp Shared library\n",
    "\n",
    "The header of the C++ shared library has this function:\n",
    "```c\n",
    "EXPORT_CODE double CONVENTION PropsSI(const char *Output, const char *Name1, double Prop1, const char *Name2, double Prop2, const char *Ref);\n",
    "```\n",
    "\n",
    "The official Windows dynamic library(64bit): `./coolpropdemo/bin/CoolProp_x64.dll`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "'J:\\\\SEU\\\\SEE\\\\PySEE\\\\home\\\\notebook\\\\demo'"
      ]
     },
     "execution_count": 19,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "%pwd"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "J:\\SEU\\SEE\\PySEE\\home\\notebook\n"
     ]
    }
   ],
   "source": [
    "%cd .. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [],
   "source": [
    "from ctypes import *\n",
    "lib = cdll.LoadLibrary(\"./coolpropdemo/bin/CoolProp_x64.dll\")\n",
    "PropsSI = lib.PropsSI\n",
    "PropsSI.argtypes = (c_char_p, c_char_p,c_double,c_char_p, c_double, c_char_p)\n",
    "PropsSI.restype = c_double\n",
    "\n",
    "def pt2h(p,t):\n",
    "    result = PropsSI(b\"H\", b\"T\",t+273.15, b\"P\",p*1.0e6, create_string_buffer(b\"Water\"))/1000\n",
    "    return  result\n",
    "\n",
    "def pt2s(p,t):\n",
    "    result = PropsSI(b\"S\", b\"T\",t+273.15, b\"P\",p*1.0e6, create_string_buffer(b\"Water\"))/1000\n",
    "    return  result"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "3400.8743470814416"
      ]
     },
     "execution_count": 22,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "pt2h(16,536)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "6.4355597642963245"
      ]
     },
     "execution_count": 23,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "pt2s(16,536)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Reference\n",
    "\n",
    "* ctypes http://docs.python.org/3/library/ctypes.html\n",
    "\n",
    "* [Numerical Recipes](http://numerical.recipes)\n",
    "\n"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.8.7"
  },
  "toc": {
   "base_numbering": 1,
   "nav_menu": {},
   "number_sections": false,
   "sideBar": true,
   "skip_h1_title": false,
   "title_cell": "Table of Contents",
   "title_sidebar": "Contents",
   "toc_cell": false,
   "toc_position": {
    "height": "calc(100% - 180px)",
    "left": "10px",
    "top": "150px",
    "width": "165px"
   },
   "toc_section_display": true,
   "toc_window_display": true
  },
  "widgets": {
   "state": {},
   "version": "1.1.2"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}
